from pulp import LpProblem, LpVariable, LpMaximize, lpSum, value, LpMinimize

import numpy as np
import matplotlib.pyplot as plt

from scipy.optimize import linprog

plt.rcParams['font.sans-serif'] = ['Simhei']  # 显示中文
plt.rcParams['axes.unicode_minus'] = False  # 显示负号


def solve_1():
    """
    图解法求解线性规划问题
    max z = 2*x1 + 2*x2
    约束条件:
        5*x1 + 2*x2 <= 15
        6*x1 + 2*x2 <= 24
        x1 + x2 <= 5
        x1 >= 0, x2 >= 0
    """
    # 生成x2的值范围 (横坐标)
    x2 = np.linspace(0, 8, 400)  # 必要

    # 计算约束边界线 (转换为x1关于x2的表达式)
    y1 = (15 - 2 * x2) / 5  # 5*x1 + 2*x2 <= 15 <=> x1 = (15 - 2*x2)/5
    y2 = (24 - 2 * x2) / 6  # 6*x1 + 2*x2 <= 24 <=> x1 = (24 - 2*x2)/6
    y3 = 5 - x2  # x1 + x2 <= 5 -> x1 = 5 - x2

    # 创建图形
    plt.figure(figsize=(7, 6))

    # 绘制约束线
    plt.plot(y1, x2, label=r'$5x_1 + 2x_2 \leq 15$')
    plt.plot(y2, x2, label=r'$6x_1 + 2x_2 \leq 24$')
    plt.plot(y3, x2, label=r'$x_1 + x_2 \leq 5$')
    plt.xlim(0, 6)
    plt.ylim(0, 6)
    plt.xlabel(r'$x_1$')
    plt.ylabel(r'$x_2$')
    # 填充可行域feasible_domain
    feasible_domain = np.minimum(np.minimum(y1, y2), y3)
    plt.fill_betweenx(x2, 0, feasible_domain,
                      where=(x2 <= 10),
                      color='grey',
                      alpha=0.5)

    # 使用线性规划求解最优解,目标函数取反，linprog求解的是最小化问题
    c = [-2, -2]  # 目标函数系数
    A = [[5, 2], [6, 2], [1, 1]]  # 不等式约束系数
    b = [15, 24, 5]  # 约束右侧常数
    bounds = [(0, None), (0, None)]  # 变量边界

    res = linprog(c, A_ub=A, b_ub=b, bounds=bounds)
    if res.success:
        # 绘制目标函数的最优解点
        plt.plot(res.x[0], res.x[1], 'ro', label='最优解')
        plt.legend()
        plt.show()

        print(f"结果: {res.fun * -1:.2f}")  # 最大化目标函数的最优值
        print(f"最优解: x1 = {res.x[0]:.2f}, x2 = {res.x[1]:.2f}")
    else:
        print("没有发现解")


def solve_2():
    """
    min z = 6x1+5x2+4x3

    x1+x2+x3 == 100
    400x1+600x2+400x3>=48000 --> -400x1-600x2-400x3 <= -48000
    800x1+200x2+400x3>=44000 --> -800x1-200x2-400x3 <= -44000
    x1,x2,x3>=0
    :return:
    """

    z = np.array([6, 5, 4])  # 确定目标函数系数

    """
        确定不等式约束(系数)，需要转换成 <= 约束
        A 为含决策变量的系数矩阵，b为右端约束
        有两个约束，依次对应
    """
    A_ub = np.array([
        [-400, -600, -400],  # -400x1-600x2-400x3 <= -48000
        [-800, -200, -400]  # -800x1-200x2-400x3 <= -44000
    ])

    b_ub = [-48000, -44000]

    # 确定等式约束条件，只有一个
    A_eq = np.array([
        [1, 1, 1]  # x1+x2+x3 == 100
    ])

    b_eq = np.array([100])

    # 决策变量不能为负的约束
    bounds = [(0, None), (0, None), (0, None)]

    # 求解
    res = linprog(z, A_ub=A_ub, b_ub=b_ub,
                  A_eq=A_eq, b_eq=b_eq, bounds=bounds,
                  method='highs')
    # 输出结果
    if res.success:
        print("各食物的最优供应量（kg）:")
        print(f"食物 P: {res.x[0]:.2f} kg")
        print(f"食物 Q: {res.x[1]:.2f} kg")
        print(f"食物 R: {res.x[2]:.2f} kg")
        print(f"最小成本: {res.fun:.2f} 元")
    else:
        print("未找到最优解")


def solve_3():
    """

    :return:
    """
    # 模型参数设置
    """
    问题描述：
        某厂生产甲乙两种饮料，每百箱甲饮料需用原料6千克、工人10名，获利10万元；每百箱乙饮料需用原料5千克、工人20名，获利9万元。
        今工厂共有原料60千克、工人150名，又由于其他条件所限甲饮料产量不超过8百箱。
        （1）问如何安排生产计划，即两种饮料各生产多少使获利最大？
        （2）若投资0.8万元可增加原料1千克，是否应作这项投资？投资多少合理？
        （3）若不允许散箱（按整百箱生产），如何安排生产计划，即两种饮料各生产多少使获利最大？
        （4）若不允许散箱（按整百箱生产），若投资0.8万元可增加原料1千克，是否应作这项投资？投资多少合理？
    """

    # 问题 1：
    """
    问题建模：
        决策变量：
            x1：甲饮料产量（单位：百箱）
            x2：乙饮料产量（单位：百箱）
        目标函数：
            max fx = 10*x1 + 9*x2
        约束条件：
            6*x1 + 5*x2 <= 60
            10*x1 + 20*x2 <= 150            
            x1, x2 >= 0，x1 <= 8

    """
    ProbLP1 = LpProblem("问题一", sense=LpMaximize)  # 定义问题 1，求最大值
    x1 = LpVariable('x1', lowBound=0, upBound=8, cat='Continuous')  # 定义 x1
    x2 = LpVariable('x2', lowBound=0, cat='Continuous')  # 定义 x2
    ProbLP1 += (10 * x1 + 9 * x2)  # 设置目标函数 f(x)
    ProbLP1 += (6 * x1 + 5 * x2 <= 60)  # 不等式约束
    ProbLP1 += (10 * x1 + 20 * x2 <= 150)  # 不等式约束
    ProbLP1.solve()
    print(ProbLP1.name)  # 输出求解状态

    for v in ProbLP1.variables():
        print(v.name, "=", v.varValue)  # 输出每个变量的最优值
    print("F1(x) =", value(ProbLP1.objective))  # 输出最优解的目标函数值

    # 问题2：
    """
        决策变量：
        x1:甲饮料产量
        x2:乙饮料产量
        x3:增加投资
        目标函数：
        max fx = 10*x1+9*x2-x3
        约束条件：
        6*x1+5*x2<=60+x3/0.8
        10*x1+20*x2<=150
        x1,x2,x3>=0,x1<=8
    """
    ProbLP2 = LpProblem("问题2", sense=LpMaximize)
    x1 = LpVariable('x1', lowBound=0, upBound=8, cat='Continuous')
    x2 = LpVariable('x2', lowBound=0, upBound=None, cat='Continuous')
    x3 = LpVariable('x3', lowBound=0, upBound=None, cat='Continuous')

    ProbLP2 += (10 * x1 + 9 * x2 - x3)
    ProbLP2 += (6 * x1 + 5 * x2 - x3 / 0.8 <= 60)
    ProbLP2 += (10 * x1 + 20 * x2 <= 150)
    ProbLP2.solve()
    print(ProbLP2.name)

    for v in ProbLP2.variables():
        print(v.name, "=", v.varValue)
    print("F2(x) = ", value(ProbLP2.objective))

    # 问题3：
    """
        在问题1的基础上进行整数规划
          问题建模：
        决策变量：
            x1：甲饮料产量（单位：百箱）
            x2：乙饮料产量（单位：百箱）
        目标函数：
            max fx = 10*x1 + 9*x2
        约束条件：
            6*x1 + 5*x2 <= 60
            10*x1 + 20*x2 <= 150            
            x1, x2 >= 0，x1 <= 8
            x1,x2为整数
    """
    ProbLP1 = LpProblem("问题一", sense=LpMaximize)  # 定义问题 1，求最大值
    x1 = LpVariable('x1', lowBound=0, upBound=8, cat='Integer')  # 定义 x1
    x2 = LpVariable('x2', lowBound=0, cat='Integer')  # 定义 x2
    ProbLP1 += (10 * x1 + 9 * x2)  # 设置目标函数 f(x)
    ProbLP1 += (6 * x1 + 5 * x2 <= 60)  # 不等式约束
    ProbLP1 += (10 * x1 + 20 * x2 <= 150)  # 不等式约束
    ProbLP1.solve()
    print(ProbLP1.name)  # 输出求解状态

    for v in ProbLP1.variables():
        print(v.name, "=", v.varValue)  # 输出每个变量的最优值
    print("F1(x) =", value(ProbLP1.objective))  # 输出最优解的目标函数值

    # 问题4：
    """
        在问题2的基础上整数规划
        决策变量：
            x1:甲饮料产量
            x2:乙饮料产量
            x3:增加投资
        目标函数：
            max fx = 10*x1+9*x2-x3
        约束条件：
            6*x1+5*x2<=60+x3/0.8
            10*x1+20*x2<=150
            x1,x2,x3>=0,x1<=8
            x1,x2为整数
    """
    ProbLP2 = LpProblem("问题2", sense=LpMaximize)
    x1 = LpVariable('x1', lowBound=0, upBound=8, cat='Integer')
    x2 = LpVariable('x2', lowBound=0, upBound=None, cat='Integer')
    x3 = LpVariable('x3', lowBound=0, upBound=None, cat='Integer')

    ProbLP2 += (10 * x1 + 9 * x2 - x3)
    ProbLP2 += (6 * x1 + 5 * x2 - x3 / 0.8 <= 60)
    ProbLP2 += (10 * x1 + 20 * x2 <= 150)
    ProbLP2.solve()
    print(ProbLP2.name)

    for v in ProbLP2.variables():
        print(v.name, "=", v.varValue)
    print("F2(x) = ", value(ProbLP2.objective))


def solve_4():
    # 创建最大化问题
    problem = LpProblem(name='Knapsack_problem', sense=LpMaximize)
    # 定义目标函数系数
    profits = [540, 200, 180, 350, 60, 150, 280, 450, 320, 120]
    # 定义约束系数
    weights = [6, 3, 4, 5, 1, 2, 3, 5, 4, 2]
    max_weight = 30
    # 定义决策变量0-1
    x = [LpVariable(name=f"x{i + 1}", cat='Binary') for i in range(10)]
    # print(x)

    # 定义目标函数
    problem += lpSum(profits[i] * x[i] for i in range(10)), 'total_profit'
    # 添加约束条件
    problem += lpSum(weights[i] * x[i] for i in range(10)) <= max_weight, 'weight_constraint'

    # 求解问题
    problem.solve()

    print('选中的物品分别是：（0表示未选，1表示选中）')
    selected_items = [int(value(x[i])) for i in range(10)]
    print(selected_items)
    print('最大利润：')
    print(value(problem.objective))


def solve_5():
    problem = LpProblem('Assignment_problem', LpMinimize)
    # 定义目标函数系数
    c = [66.8, 75.6, 87, 58.6, 57.2, 66, 66.4, 53, 78, 67.8, 84.6,
         59.4, 70, 74.2, 69.6, 57.2, 67.4, 71, 83.8, 62.4]
    # 定义决策变量
    x = [LpVariable(f"x{i + 1}", cat='Binary') for i in range(20)]
    # 定义目标函数（最小化总时间）
    print(c[0])

    problem += lpSum(c[i] * x[i] for i in range(20)), 'total_time'

    # 添加不等式约束，每个人只能选择四种泳姿之一
    A = [
        [1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
        [0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
        [0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0],
        [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0],
        [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1],
    ]
    b = [1, 1, 1, 1, 1]
    for i in range(len(A)):
        problem += lpSum(A[i][j] * x[j] for j in range(20)) <= b[i], f'problem_constraint{i + 1}'

    # 添加等式约束
    Aeq = [
        [1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0],
        [0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0],
        [0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0],
        [0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1],
    ]
    beq = [1, 1, 1, 1]
    for i in range(len(Aeq)):
        problem += lpSum(Aeq[i][j] * x[j] for j in range(20)) == beq[i], f'style_constraint{i + 1}'

    # 求解
    problem.solve()
    print('分配结果：（0表示未选，1表示选钟）')
    assignments = [int(value(x[i])) for i in range(20)]
    print(assignments)
    print("最小时间：")
    print(value(problem.objective))

    assignment_matrix = np.array(assignments).reshape(5, 4)
    print(assignment_matrix)


def test():
    z = np.array([-20, -30, -45])  # 确定目标函数系数
    A_ub = np.array([
        [4, 8, 15],  # -400x1-600x2-400x3 <= -48000
        [1, 1, 1]  # -800x1-200x2-400x3 <= -44000
    ])

    b_ub = [100, 20]

    # 决策变量不能为负的约束
    bounds = [(0, None), (0, None), (0, None)]

    # 求解
    res = linprog(z, A_ub=A_ub, b_ub=b_ub,
                  bounds=bounds,
                  method='highs')
    # 输出结果
    print(-res.fun)


if __name__ == '__main__':
    solve_5()